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,f SECTION I 

INTRODUCTION 

Existing snowmelt runoff models are used to make single— valued 
estimates for a specified decision period. Such models are based on 
the expected value theorem for statistical decisionmaking.. A single 
valued estimate fails to provide a measure of the accuracy of an 
estimate. While a computed standard error may be used to indicate the 
accuracy, existing decisionmaking techniques that are based on the 
expected value theorem do not provide a mechanism for altering the 
decision based on the computed standard error. That is, existing 
decisionmaking techniques are not sensitive to the accuracy of the 
estimated value of snowmelt runoff. Models based on the expected 
value theorem also fail to provide for inaccuracies in the input data, 
especially the empirical coefficients that are specific to a particular 
watershed. Because these coefficients cannot be estimated with 
certainty for the specified decision period and may vary with time, 
error in the estimated values of the coefficients will cause an error 
in the output function, i.e., the runoff forecast. Incorrect 
decisions may result from inaccuracies in the output function. 

The snowmelt-runoff model is a means of forecasting snowmelt run- 
off volumes using a small number of variables. The model first 
estimates the volume of snowmelt for each day of the forecast period 
and then estimates the amount of melt water that will leave the water- 
shed on each day by using a recession curve. This model has been 
calibrated on a number of mountainous watersheds. 
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The depth of snowmelt for any day is estimated in the model by 
multiplying the number of degrees above freezing by a degree-day 
factor to obtain a depth measurement. This depth is converted into 
volume by multiplying the depth estimate by the area of the snowpack. 

The volume of precipitation for the day is then added m to get an 
estimate of the total amount of water available for runoff. This water 
is assumed to leave the watershed at a decaying rate. The model for 
predicting the snowmelt runoff on day n+1 is given by: 

Q_ , 1 = C fa (T +AT )S +P 1 (1-k n ) +Q k .. (1-1) 

n+1 n|_ n n n n nj 86400 n+1 n n+1 

in which Q is the average daily discharge (m^/S) ; C is the runoff 
coefficient; a is the degree-day factor (cm* °C - *d - ^) indicating, the 
snowmelt depth resulting from 1 degree-day; T is the number of degree- 
days (°C*d); AT = the adjustment by temperature lapse rate necessary 
because of the altitude difference between the temperature station and. 
the average hypsometric elevation of the basin or zone; S is the ratio 
of the snow-covered area to the total area; P is the precipitation 

2 

contributing to runoff (cm) ? A is the area of the basin or zone (m ) ; 

2 —*1 

0.01/86400 is a conversion from CM*m *d~ to m /S; k is the recession 
coefficient indicating the decline of discharge in a period without 
snowmelt or rainfall; and n is the notation indicating the day for which 
the- value is given. The important parameters of the model that must be 
calibrated for any specific basin are the degree-day factor, the recession 
constant, and a runoff coefficient that is included to account for 
losses due to infiltration and evapotranspiration. Additionally, the 
value of AT must be considered as a variable because its value is not 
known exactly; it is based on an extrapolation and a mean zonal 
hypsometric elevation. 


Project Objective 


The goal of this project is to develop a stochastic form of 
the snowmelt runoff model that can be used for probabilistic decision- 
making. The use of probabilistic streamflow predictions instead of 
single-valued deterministic predictions should lead to greater 
accuracy in decisions. While the accuracy of the output function 
is important in decisionmaking, it is also important to understand 
the relative importance of the coefficients. Therefore, a sensitivity 
analysis will be made for each of the coefficients. 
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SECTION II 
LITERATURE REVIEW 

Many different models have been developed and used for 
predicting snowmelt runoff (Leaf, 1977; Baker and Carder, 

1977; Zuzel and Cox, 1978). These models vary considerably 
in complexity; the simplest models are based solely on sta- 
tistical techniques, while the most complex methods attempt 
to model the individual processes involved in the melting of 
a snowpack. Some models are designed to predict streamflow 
for any given day or series of days, (Leaf, 1977; Martinec, 

1975; Tangborn, 1977) while other models give only seasonal 
predictions (Zuzel and Cox, 1978). Generally, snowmelt models 
may be categorized on the basis of complexity and length of 
forecast period.. 

Empirical models are based on statistical correlations 
between predictor variables and the criterion variable, volume 
of snowmelt runoff. This type of model is most often used for 
seasonal predictions. Snow water equivalent measurements, prerr 
vious runoff volumes, and precipitation totals are the most, 
common predictor variables (SCS, 1970; USACE, 1956). Theory 
is not very important in formulating empirical models; the ob- 
jective is to explain as much of the variation in the criterion 
values as possible using whatever data are available. It is 
quite common for these models to include two predictors express- 
ed in different units, such as snow water equivalent (in inches) 



I 


5 


and previous winter runoff (in volumetric units). 

Water balance models are more conceptual than the simple 
empirical models. The water balance is an accounting of all 
the water entering and leaving the basin. The volume of water 
stored in the snowpack is estimated from precipitation or water 
equivalent data; allowances are made for losses due to evapora- 
tion t groundwater storage, and transpiration; the remaining volume 
is the seasonal snowmelt runoff prediction (Zuzel and Cox, 1978). 
Loss rates may be estimated either empirically or conceptually, 
as may the snowpack storage. Most water balance models are some- 
what empirical. 

Short-term runoff predictions usually require models of great- 
er complexity than the models used for seasonal runoff. Not- only 
must the. total volume of water stored in the snowpack be estimat- 
ed, but also the proportion of that volume that will melt and 
leave the watershed as streamflow in a given time period must 
be estimated. The amount of water generated by melting snow 
is a function of the energy available for this purpose. There- 
fore, the most complex snowmelt models are generally based on an 
energy balance (Zuzel and Cox, 1978). 

Energy balance procedures attempt to model the physical pro- 
cesses involved in snowmelt runoff. The amount of available energy 
is commonly estimated by the air temperature, although some models 
include such factors as incoming solar radiation, cloud cover, 
albedo, and net long-w.ave radiation (Anderson, 1976). These 
models often require that the watershed be subdivided into small, 
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homogeneous areas so that the available energy for each location 
can be estimated more accurately (Leaf, 1977). Since snowmelt 
models are generally used in mountainous areas, slope and aspect 
can result in large differences in incident energy from one area 
to another. Evaporation, transpiration and groundwater losses 
are also estimated conceptually in some energy budget models 
(Leaf, 1977). 

Model Selection 

To test the study objectives, models having significant dif- 
ferences in important characteristics had to be selected. Cri- 
teria for model selection include the frequency of current usage, 
input data requirements and whether or not these data are typi- 
cally available, the degree of model complexity, and the length 
of forecast period. Additionally, because snow covered area 
(SCA) is more readily available than in previous decades, models 
that either included SCA or were capable of being modified to 
include it were given more consideration. 

Three models were selected for comparison, with several 
methods of evaluation for each model. The model types studied 
were the regression model, the Tangborn model, and the Martinec 
model. 

The Regression Models 

The most common form of empirical model is the linear re- 
gression. These models are widely used for snowmelt runoff pre- 
dictions in the western U.S. (USACE, 1956; SCS, 1970). They are 
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easily calibrated and can use many different hydrologic vari- 
ables as predictor variables. These models are used for mak- 
ing seasonal runoff forecasts, but due to the empirical nature 
of the method, they may also be used to give predictions for 
shorter time periods. 

Linear regression models are based on the assumption that 
there is a linear relationship between the predictor variables 
and the criterion variable. This assumption implies that as 
the value of the predictor variable increases, the value of 
the criterion variable changes at a constant rate. The equa- 
tion that relates the value of the criterion to the value of 
the predictor is of the form: 

Y - a + bX (2-1) 

in v/hich Y is the criterion variable, X is the predictor vari- 
able, and a and b are the regression coefficients (Miller and 
Freund, 1977). 

Many hydrologic variables have approximately linear rela- 
tionships with the volume of snowmelt runoff. A few of these 
variables are snow water equivalent, winter precipitation, and 
snow covered area. The linearity of the relationships is due 
to the fact that these variables are indicators of the volume 
of water stored in the snowpack. Because the relationships 
between these predictor variables and the volume of runoff 
are only approximately linear, many different lines may be 
drawn which appear to fit the data. Some of the lines pass 
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through a number of the data points , but due to deviations 
from linearity, a straight line that will pass through all 
of the data points can not be drawn. 

The method of selecting the best regression line for a 
set of data points is based on minimizing the sum of squares 
of the errors. For each observed value of the predictor, two 
values of the criterion variable appear; the first is the cor- 
responding observed value and the second is the value predict- 
ed by the regression equation. The difference between these 
two values is termed the error of prediction. The regression 
line is defined as the line that results in the minimum value 
of the sum of the squares of the errors. The coefficients of 
the regression line can be derived using the, equations : 

b _ IXY - (ZXEY)/n 
EX 2 - (LX) 2 /2 

and 

a = (ZY)/n - b(EX)/n 

in which X and Y are the predictor and criterion variables, 
respectively, and n is the number of observations (Hays, 1965). 

By using these equations, the line of best fit can be determined. 

In natural systems the value of the criterion variable is 
often a function of more than one predictor. The relationships 
between the criterion variable and the predictors may be assumed 
to be linear, resulting in a prediction equation of the form: 


(2-2) 


(2-3) 
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Y = a + I b.X. (2-4) 

i=l 1 1 

t h. 

in which Y is the criterion variable, is the i L predictor 
variable, and a and b^ are the regression coefficients. Models 
of this type are called multiple linear regressions. The re- 
gression coefficients are unique and may be calculated from 
equations similar to Eqs. 2-2 and 2-3. In many cases, the in- 
clusion of more than one predictor variable results in a more 
accurate model (Davis, 1973)* 

The Tangborn Model 

The Tangborn equation is a water balance model (Tangborn 
and Rasmussen, 1976). The structure of the model was establish- 
ed conceptually, but calibration is accomplished using regression 
methods. The model may be used for any length of forecast period 
from one day to the entire snowmelt season. The only data re- 
quired are daily precipitation and runoff values , although daily 
temperature may be included for short forecast periods. 

The basic form of the model is: 

R* = a P w + b - R w (2-5) 

in which R* is the predicted runoff volume, P y is the total 
depth of precipitation observed during the preceding winter, 

R is the winter runoff, and a and b are regression coefficients. 
The structure of the model is based on the assumption that the' 
volume of water stored on the watershed is equal to the amount 
of winter precipitation minus the winter runoff. The regression 
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coefficients represent losses and modifications such as trans- 

fit 

piration, groundwater storage, and evaporation. 

An important feature of the Tangborn model is the test 
season modification. In using this method, a short test season 
prediction model with the structure of Eq. 2-5 is developed. 

At the end of the test season, the error of the test season 
prediction is evaluated and used to modify the prediction for 
the forecast season. The form of the forecast model becomes: 

R** = R* - ce. = a(P +P. ) + b - (R +R. ) - ce. (2-6) 

S S 0 Wu Wt u 

in which R is the revised runoff prediction; R is the original 
prediction; P w and P^ are the winter and test season precipita- 
tion, respectively; R w and R t are the winter and test season run- 
off volumes, respectively; a, b, and c are coefficients; and e fc 
is the error of the test season prediction. The reasoning be- 
hind this modification is that the test season error is a result 
of the inaccuracy of estimating basin storage by subtracting 
winter runoff from winter precipitation. Because the forecast 
season prediction is based on the same estimate, the test season 
error should be related to the prediction season error. 

Figure 1 shows the relationship of the various seasons. In 
order to use the test season approach, data from the present and 
a number of previous years are compiled. For- each year, precipi- 
tation and runoff totals are computed for the winter and test 
seasons; runoff totals are also computed for the prediction 
season of each year, except for the current year (the value for 
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the current year is not yet known). Note that the prediction 
date, April 1, is at the end of the test season; therefore, 
observed values of runoff and precipitation during the test 
season are available for the current year. Once all the data 
has been obtained, the observed test season runoff volumes are 
regressed onto the winter precipitation values, resulting in a 
calibrated equation of the form: 

< = a p „ + b - (2 -?> 

♦ 

in which is the predicted test season runoff. The test season 
error in each year is then computed by the equation: 

«* = *£- R t (2 - 8) 

Next, a model for estimating the prediction season runoff is form 
ed by regressing the prediction season runoff on the sum of the 
winter and test season precipitation depths for each of the pre- 
vious years: 


R* = a(P w +P t ) + b - (R w +R t ) (2-9) 

The errors are then calculated in a manner similar to that used 
for the test season: 

= R s - R s (2 - 10) 

in which e g is the prediction season error. The coefficient of 
the test season error, c in Eq. (2-6), can then be determined. 



The coefficient is computed using the test season and predic- 
tion season errors from previous years, according to the equa 
tion: 
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c 


z:(e t* e s ) 

I(e t ) 2 


( 2 - 11 ) 


The original runoff season prediction, R ,. which was calcu- 
lated for the current year in Eq. (2-9)? is adjusted by the pro- 
duct of c and the current year test season error; the final pre- 
diction is: 


C = H s - ce t = a( v p t> + b - <w - ce t (2 - 12) 

in which R g is the final prediction. 

When using the Tangborn model for prediction periods of a 
few days, accuracy may be increased by including temperature in 
the model (Tangborn, 1978). Tangborn suggested the following 
composite temperature variable, A^: 

A t = a¥ + (l-a)AT (2-13) 

in which T is the daily mean temperature , AT is the daily range of 
temperature, and a is a coefficient. The daily mean temperature 
is computed from the observed maximum and minimum temperatures 
for the day; the range of temperature is the difference between 
the maximum and minimum observed values. The reasoning behind 
this equation is that the average daily mean temperature is an 


estimator of the amount of convective energy available for 
melting snow, and that the difference between maximum and mini- 
mum temperatures can be used to estimate the amount of radiant 
energy available for this purpose. Large differences between 
the daily maximum and minimum are indicative of clear skies, 
while a small daily range of temperature indicates cloud cover 
and, therefore, less radiant energy. The relative importance 
of the two components (radiative and convective) is controlled 
by the coefficient a. When the temperature tern is included 
in the Tangborn model, the equation becomes: 

jfc jfc jfc jjc 

R = a(P w +P t ) + b - (R w +R t ) - ce t - e s (2-14) 

s 

in which e* is. the prediction season error estimated from the 

20C $ 

temperature function, and R is the revised runoff prediction. 
The value of the prediction season error is estimated from the 
temperature function A t using the equation: 

e* = dA t + e (2-15) 

in which d and e are coefficients determined by regression. 
Tangborn reports a minimum reduction in standard error of 
estimate of nine percent due to inclusion of this tempera- 
ture term (Tangborn, 1978). 

The Martinec Model 

The Martinec model is conceptually derived and may be 
used for prediction periods of one day or longer (Martinec ,1975) ♦ 
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The amount of energy available for anowmelt runoff is estimated 
by a daily temperature index. Data requirements include daily 
temperature, precipitation, and snow covered area. The form of 
the model is: 

Q* = c (dTSCA+P) A(l-K) + KQ„ , (2-16) 

n n-i 

in which Q* is the predicted volume of runoff for day n, c is a 
dimensionless runoff coefficient, d is a degree-day factor, T is 
the value of the daily temperature index on day n, A is the total 
area of the watershed, SCA is the percentage of the area that is 
covered by snow on day n, K is a dimensionless recession coeffi- 
cient, and Q n-1 is the volume of runoff observed on the previous 
day. The value of the daily temperature index is computed using 
hourly data if available; otherwise, the daily maximum and minimum 
temperatures are used. The daily index is a measure of the average 
number of degrees above freezing for the temperature on that day. 

The values are expressed in degree-days Celsius. 

The first term of Eq. (2-16) represents the amount of water 
that is generated by precipitation and melting snow on day n and 
that is expected to leave the watershed on that day. The value of 
the degree-day factor, d, is expressed in inches of water per degree 
Ceisium; therefore, when the temperature index is multiplied by this 
factor, an estimate of the depth of water generated by snowmelt is 
obtained. This depth is multiplied by the total area of the water- 
shed, A, and by the percentage of the total area that is covered by 
snow (SCA) to get an estimate of the volume of water produced by 
melting snow on day n. The precipitation, P, is assumed to be a 
constant depth over the entire watershed; therefore, the product 
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of P and A is an estimate of the volume of rainfall on day n. 

The siim of the volume of melted snow and the volume of precipi- 
tation is referred to as the generated runoff. 

Not all of the generated runoff leaves the watershed on 
the day of generation. Some is lost to groundwater storage 
and evapotranspiration; this proportion is represented’ by c, 
the runoff coefficient. Furthermore, on large watersheds the 
outlet of the basin is quite a distance from the source of much 
of the generated melt; therefore, much of the water is in transit 

f 

to the outlet for several days. The proportion of water that 
does not reach the outlet on the day that it is generated is 
represented by K, the recession coefficient. Thus, only the 
proportion (1-K) of the runoff generated on day n actually reach- 
es the outlet on day n. 

The second term in the equation, K*Q n ^ , is called the re- 
cession term. It represents the amount of water generated on pre- 
vious days that is expected to appear as runoff on day n. Be- 
cause K is nearly equal to 1 on large watersheds, this recession 
term is often considerably larger than the generated runoff term. 
The Snowmelt Runoff Model 

The snowmelt -runoff model is an advanced form of the Martinec 
model. Both models can be used to forecast daily streamflow in 
basins where snowmelt is a primary contribution to the total 
runoff. The snowmelt-runoff model is designed to be used with 
remotely sensed estimates of snow covered area, as well as tempera- 
ture and precipitation data. The User's Manual (Martinec, et al., 
1983) provides information on the computer program, the required 
input, and the application of the model. 
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SECTION III 

SENSITIVITY ANALYSIS OF THE SNOWMELT-RUNOFF MODEL 
Hydrologic modeling is a procedure in which one or more of the 
phases of the hydrologic cycle are represented by a simplified system. 
Although the physical principles underlying the specific hydrologic 
processes should be considered in formulating a model structure, the 
final design invariably is only an approximation of the processes being 
modeled. In spite of the approximations involved, models can often 
provide insight into those parts of the physical process in which 
knowledge of the underlying principles is deficient. 

In addition to being familiar with the physical principles being 

s 

modeled, the designer must also be familiar with modeling tools that 
aid in formulating, calibrating, and verifying conceptual representa- 
tions of the unknown parts of the physical processes. Sensitivity 
analysis is a modeling tool that, if properly used, can provide a model 
designer with a better understanding of the correspondence between the 
model and the physical processes being modeled. Sensitivity of model 
components and parameters is potentially useful in the formation, 
calibration, and verification of a hydrologic model. However, in the 
past the use of sensitivity has been limited to the determination of an 
optimal set of model parameters and identifying the effect of variability 
in a parameter on the optimal solution. Parameter optimization tech- 
niques categorized as gradient climbing procedures use estimates of 
parameter sensitivity to derive an optimal set of model parameters. 

In post-optimization analysis, sensitivity has frequently been used as 
a means of ranking the variables in order of relative importance. 
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MATHEMATICAL FOUNDATIONS OF SENSITIVITY 


Definition 

Sensitivity is the rate of change in one factor with respect to 
change in another factor. Although such a definition is vague in 
terms of the factors involved it nevertheless implies a quotient of 
two differentials. Stressing the nebulosity of the definition is 
important because, in practice, the sensitivity of model parameters 
is rarely recognized as a special case of the concept of sensitivity. 
The failure to recognize the generality of sensitivity has been 
partially responsible for the limited use of sensitivity as a tool for 
the design and analysis of hydrologic models. 

The Sensitivity Equation 

The general definition of sensitivity can be expressed in 
mathematical form by considering a Taylor series expansion of the 
explicit function: 


0 = f(F r F 2 ,F ) 


(III-l) 


The change in factor 0 resulting from change in a factor F^ is given, 
by: 

30 „ , 3 0 „ 

(III-2) 


30 3^0 ? 

f< V* F i- — !K 

X Ox. 


1 

in which is the value of 0 at some specified level of each F^. 

If the nonlinear terms are small in comparison with the linear terms, 
Eq. III-2 reduces to: 


3 °n 

f ( F . +A F . , F . , = O rt + 

11 j|j^i 0 3F. 


•AF. 

x 


( 1 1 1 — 3 ) 
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Thus: 

AO q = f(F. + 4F., F j|j*i>- O 0 =(aO 0 /3 Fi> SF i 
Equation III-4 is referred to herein as the linearized sensitivity 
equation. It measures thechange in factor 0 that results from change 
in factor F^. The linearized sensitivity equation can be extended to 
the case where more than one parameter is changed simultaneously. The 
general definition of sensitivity is derived from Eq. III-l and III-4: 


S=30_/3F . =[f(F.+AF , F . i , ,. ) -f (F, ,F. F ) ]/AF. 


j| j 5^i 


1 ' 2 ‘ 


n 


(HI-5) 


O' x X x 

Computational Methods 

The general definition of sensitivity which is expressed in 
mathematical form by Eq. III-5 suggests two methods of computation. 

The left-hand side of Eq. III-5 suggests that the sensitivity of O 
to changes in factor F^ can be estimated by differentiating the explicit 
relationship of Eq. 1 with respect to factor F^ : 


S =30 /3F. 
o x 


(III-6) 


Analytical differentiation has not been used extensively for analyzing 
hydrologic models because the mathematical framework of sensitivity 
has not been sufficiently developed. It will be used even less, 
frequently as hydrologic models become more complex. 

The method of factor perturbation, which is the second computa- 
tional method suggested by Eq. III-5, is the more commonly used method 
in hydrologic analysis. The right-hand side of Eq. III-5 indicates 
that the sensitivity of O to change in F can be derxved by incrementing 
Fj and computing the resulting change in the solution O. The sensi- 
tivity is the ratio of the two differentials and can be expressed in 


finite difference form: 
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S=AO 0 /&F i =[f (F i +AF i/ F j | jV± )-f (F r F 2 Fj ]/ A F.. (III-7) 

However, use of the method of parameter perturbation is often imprac- 
tical for a complete sensitivity analysis of multiparameter systems 
because of the extensive computational effort required for complex 
models. It is used when analytical analysis of the sensitivity is 
not possible. 

Parametric and Component Sensitivity 

A simplified system or a component of a more complex system is 
described by three functions: the input function, the output function 

ft 

and the system response function. The response function is the 
function that transforms the input function into the output function 
and is often defined by a distribution function which depends on one 
or more parameters. In thepast, sensitivity analyses of models have 
been limited to measuring the effect of parametric variations on the 
output. Such analyses focus on the output and response functions. 

I 

Using the form of Eq. III-5 parametric sensitivity can be mathematically 
expressed as: 


30 

3p i 


f (Fj+AP^Pj i j^j) ~ £ (Pj/P 2 P n > 

AP. 

l 


(III-8) 


where 0 represents the output function and is the parameter under 
consideration. 

Unfortunately, the general concept of sensitivity has been over- 
shadowed by parametric sensitivity. As models have become more complex 
the derivation of parametric sensitivity estimates have become 
increasingly more difficult and, most often, impossible, to compute. 
However, by considering the input and output functions the general 
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definition of sensitivity, Eq. III-5, can be used to derive another 
form of sensitivity. Component sensitivity measures the effect of 
variation in the input function I on the output function: 


30 _ AO 
C 31 AI 


(III-9) 


Combining component and parameter sensitivity functions makes 
it feasible to estimate the sensitivity of parameters of complex models. 
For example, in the simplified two component model of Fig. III-l, the 
sensitivity of Y to variation in and the sensitivity of Z to 
variation in P 2 are readily computed using sensitivity as defined by 


Eq. IIX-6 : 

= 3Y/3P 1 and S 2 = 3Z/3P 2 (HI-10) 

However, the sensitivity of the output from component 2 to change in 
the parameter of component 1 cannot always be estimated directly from 
the differential 3Z/3P^. In such cases, the component sensitivity 
function of component 2 can be used with the parametric sensitivity 
function to estimate the sensitivity of Z to change in P^. 

Specifically, the sensitivity of 3Z/3P.J, equals the product of the 
•component sensitivity function 3Z/3Y and the parametric sensitivity 
function 3Y/3P^: 

az/aPj = (3Z/3Y) • Oy/sp^ (iii-ii) 


Whereas the differentials 3Z/3Y and SY/SP^ are often easily derived, 
an explicit sensitivity function 3Z/3P^ can be computed only for very 
simple models. When a solution cannot be obtained analytically, then 
the numerical method of Eq. III-7 must be used. 




FIGURE III-l. Sensitivity Analysis of Two-Component Model 
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Absolute and Relative Sensitivity 

Sensitivity can be expressed in two forms: absolute and relative. 

The form in which sensitivity values are presented depends on the 
intended use. Sensitivity values computed with the definition of 
Eq. III-5 are in absolute form. Such a definition is inappropriate 
for the comparison of sensitivity values because values computed 
using Eq. III-5 are not invariant to the magnitude of either factor 
0 or F^. Dividing the numerator of Eq. III-5 by 0 Q and the denominator 
by provides an estimate of the relative change in 0 with respect 
to a relative change in F^ : 


30/0 q 

R S = 3F./F. 

^ r -l 


30 

3F . 
X 


F. 

1 


(III-12) 


I X x o 

Parametric Sensitivity in the Optimization Process 

Parametric sensitivity is a vital part of all calibration 
strategies. It is used in analytical optimization when the derivatives 
of the objective function with respect to the unknowns are taken. 

For example, the principle of least squares computes the derivative 
of the sum of the squares of the errors with respect to each of the 
regression coefficients. In Lagrangian optimization, it is also 
necessary to take derivatives with respect to coefficients in the 
constraint (s) . 

In addition to analytical optimization, sensitivity analysis 
plays a central role in numerical optimization. A set of parametric 
sensitivity coefficients defined by Eq. III-5 represents the gradient 
of the output function. The parametric response surface defined by. 
the output function of multiparameter hydrologic models invariably 
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contains a considerable number of stationary points (points at which 
the gradient is zero). Thus, the gradient techniques which are useful 
for locating a stationary point in a localized region are mathematically 
insufficient for locating a global optimum. Selecting the global 
optimum from among the many stationary points requires a more complete 
optimization strategy. However, parametric sensitivity will be a 
vital component of any strategy selected for numerical optimization. 

Systematic optimization techniques based on the gradient of the 
output function in the parametric hyperspace are rarely used with 
complex simulation models. Use of such techniques are limited by the 
inability to efficiently estimate the sensitivity of many parameters. 

For models involving as few as five parameters that require evaluation, 
the method of parameter perturbation requires a considerable amount 
of computer time to evaluate the sensitivity of each parameter. 

Component sensitivity may, in the future, demonstrate the efficiency 
of the direct method of differentiation (Eq. III-6) for estimating 
the optimal values of unknown parameters. 

Stability of the Optimum Solution 

A sensitivity plot is a graphical comparison of the percent change 
in output and the percent change in a parameter value. The change in 
the value of an objective function is often used to represent the 
change in output. The sensitivity plot can be used to examine the 
stability of a parameter of the optimum solution. Derivation of the 
sensitivity plot has traditionally involved an iterative procedure 
in which the percentage change in the value of the objective function 
or output is computed for different percentage changes in a parameter 
value. For multi-parameter models the required computer time Is often 



excessive. The direct method of differentiation (Eq. III-6) is an 
alternate means of examining optimum solution stability having the 
advantage of requiring significantly less computer time. 

The sensitivity function approach is strictly valid only when 
the nonlinear terms of Eq. III-2 are insignificant. However, even 
for large changes in a parameter value the sensitivity function 
approach provides a reliable qualitative indication of the stability 
of the optimum solution. Furthermore, the sensitivity function 
approach requires considerably less computational effort than the 
derivation of a sensitivity plot. 

Sensitivity of Initial Value Estimate 

An estimate of the initial state of the system is often necessary. 
For example, simulation models involving storage of water in snow 
pack or soil zones require estimates of the initial water content of 
each storage component. Also, models for estimating watershed reten-- - 
tion rates often require an estimate of the initial rate. If such 
models are to be used for design purposes, then it is desirable to 
estimate the effect of error in initial value estimates- on the 
design variables. A sensitivity analysis of a proposed initial estimate 
can be used to estimate the effect of error in initial estimates on 
the computed output. 

The method of parameter perturbation is usually used to measure 
the effect of error in initial value estimates. However, the 
differential approach to estimating sensitivity requires considerably 
less computational effort than the perturbation technique, especially 
for complex models. 
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Sensitivity and Data Error Analysis 

Data used for model verification or for estimating values of model 
parameters invariably contain error. The magnitude of distribution 
of error often cannot be evaluated. For example, air speed is often 
included in evaporation models as a measure of air instability. But 
whether or not air speed measurements quantitatively measure the 
effect of air instability on evaporation rates is difficult to assess. 
However, the effects of data error from other sources can be 
quantitatively evaluated. For example, instrument specifications- 
supplied by manufacturers can be used to estimate the potential error 
in an observed measurement due to the inaccuracy of the recording 
instrument. Reproducibility is a measure of how well different 
testers can estimate the value of a property with the same input 
information and procedure. The lack of reproducibility is a measure 
of random error. Sensitivity and error analyses are a means of 
examining the effect of the lack of reproducibility on a hydrologic 
design variable. 

An error analysis can be made using the linearized sensitivity 
equation of Eq. III-4. In this case, the derivative 30 q/ 9F^ is the 
sensitivity of O with respect to F^ , af^ is the error in F^, and 
AOq is the error in Oq that results from AF^. One convenient method 
is to use the standard error of a factor as the error variation AF^. 
Thus, Eq. I I I- 4 becomes: 

A0 n = (30 o /3F. )S . (HI-13) 

u u i ei 

in which S e ^ is the standard error of F^. If Eq. III-13 is computed 
for each factor, then the values of AOq for each factor can be computed 
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and compared. For example, if S ^ were the standard errors of the 
measurement error, then the computed value of for each would 

indicate the relative importance of the measurement error in the 
prediction of the output O. 

Limitations of Sensitivity Analysis 

Sensitivity analysis was shown to be a useful tool for all 
phases of modeling (formulation, calibration and verification) , as 
well as part of an error analysis for decisionmaking. However, it 
has some limitations. First, sensitivity analysis is usually applied 
using the linear sensitivity equation (Eq. III-4) ; however, the 
linear form is valid only over a limited range of the variable in. 
the denominator. Since most hydrologic models are nonlinear, including 
the snowmelt-runoff model, the sensitivity changes as the value of 
each parameter changes. In most applications, sensitivity coefficients 
are usually computed using the means of the variables as a base point. 
However, the sensitivity at the extreme values of the physical 
conditions may also be of primary interest. 

The univariate nature of sensitivity is a second limitation. 

In general, sensitivity functions are derived while hodling the values 
of the other variables constant at some base point value. Similarly, 
for the error analysis the reproducibility errors represent the 
"independent" effects of the input variables. That is, the repro- 
ducibility errors assume error only in the single variable of interest. 
An underlying assumption is that there is no interaction between the 
variables. This is usually not true. Also, we cannot assume that 
any user on one design will make an error on only one variable. 
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However, because of the interaction between variables one cannot just 
simply add the reproducibility errors. For example, even if a user 
made an error of one standard error in each of the variables, the 
net affect would not be the sum of the reproducibility errors. 

Thus, the univariate nature of sensitivity is a limiting factor for 
large errors or deviations in the variables. 

The third, and probably the most important, limitation of 
sensitivity analysis is that it provides only a single-valued 
indication of the effect on the criterion or dependent variable. 

Ideally, one would like to have some idea of the distribution of the 
design variable. One could approximate this by using the linearized 
sensitivity equation with the distribution of the error or variation 
of the indpendent or input variable, but this approach is limited in 
usefulness because of the first two limitations, i.e., linearity and 
univariate. A method that circumvents the limitations of sensitivity 
analysis would be an improvement for many analyses of design methods. 

This does not imply that senitivity is not of value, only that advanced 
forms of sensitivity analysis are needed. The linear, univariate 
form of sensitivity analysis will still be of value for analyses 
described previously especially for analytical and numerical optimization 

APPLICATION OF SENSITIVITY ANALYSIS TO THE SNOWMELT- RUNOFF MODEL 

The general structure of the snowmelt-runoff model is given by: 


Q .,=C fa (T +AT )S +P ( 1-k _) + Q k , (III-14) 

n+1 n n n n n n 86400 n+1 n n+1 

in which Q is the average daily discharge (m^/S) ; C is the runoff 
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coefficient expressing the losses as a ratio (runoff/precipitation) ; 
a is the degree-day factor (CM/°Cd) indicating the snowmelt depth 
resulting from 1-degree-day; T is the number of degree-days (°C*d); 

AT is the adjustment by temperature lapse rate necessary because of 
the altitude difference between the temperature station and the 
average hypsometric elevation of the basin or zone; .S is the ratio 
of the snow covered area to the total area; P is the precipitation 

2 

contributing to runoff (CM) ; A is the area of the basin or zone (M ) ; 

k is the recession coefficient indicating the decline of discharge 

in a period without snowmelt or rainfall (i.e., k = Q ,,/Q ); n is the 

n+i n 

sequence of days during the discharge computation period; and 
0.01/86400 is a conversion from CM*M 2 /d to M'Vs. Eq. III-14 assumes 
that there is a single watershed unit; for watersheds with nonhomogeneity 
of runoff, rainfall, or temperature characteristics, the watershed 
should be subdivided into homogeneous subareas. When the watershed is 
subdivided, Eq. III-14 must be modified to reflect the subdivision. 

This requires repeating the first term of the right-hand side of 
Eq. III-14 for. each subarea. 

The model of Eq. III-14 contains eight quantities. The drainage 
area is constant for a watershed or subwatershed. Time-varying input 
variables include T, S, and P; these will vary daily. The remaining 
four quantities (a, C, AT, and k) may vary on a daily basis; however, 
they are not directly measureable, and so without additional input, 
they would probably be considered as constants for a given snowmelt 
season. The assumption that these quantities (i.e., c, a, k, and AT) 
are constant will introduce an error, which can be assessed using a 
sensitivity analysis. Rewriting Eq. III-14 yields: 
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°„ + l =lC a S n <T n +AT)+a, n 1 (1 - k,+Q n k (1II - 15) 

The subscripts for C, a, AT, and k were dropped to indicate that they 
are assumed to be constant. The derivatives of Eq. III-15 with 
respect to the four constants are: 


3Q 


n+1 


3C 


, 0. OlA 
' 86400 


) (1-k) [a 


S (T+AT) +p 1 
n n 


(III-16) 


3 Q 
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(III-18) 


3Q 


n+1 


3k 


= Q 


n 


, 0 .01A 
'86400 


) [C a S (T+AT) +CP ] 
1 n n 1 


(III-19) 


An error analysis can be performed by multiplying the derivatives of 

Eqs. III-16 to III-19 by the error of the corresponding variable 

(see Eq. III-4) . The error in Q n+ ^' which can he denoted as AQ n+ ^, 

is a function of the daily values of S, T, and P, and for computing 

the error with respect to k it would also be a function of Q n - 

To illustrate the use of Eqs. III-16 to III-19, consider the 

following values for an 8.9 sq. km. subwatershed: C=0.95, a=0.45, 

k=0. 87, S=0. 8, T=1 .15, AT=0 .65, P=0.21, and Q =0.453. Eq. III-14 

n 

yields a predicted runoff rate of 0.503 M^/S. Eqs. III-16 to III-19 
yield the following derivatives: 3Q n+ ^/3C=0 . 297, 3Q n+ ^/3a=0. 183, 

3Q n+ i/3 AT=0 . 0458, and 3Q n+ ^/3k=-0 . 3866 . These values cannot be 
compared to measure the relative importance of the four quantities 
(C, a, AT, and k) . 
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In order to use the derivatives as part of an error analysis, 
values of the change in C, a, AT, and k must be estimated. Based on 
the data for the Dischma watershed given in the Snowmelt-Runoff Model 


User's Manual the following errors were estimated: AC=0.025, 


Aa=0.05, A(AT)=0.01, and Ak=0.15. The errors in Q that would 

n+1 

result from these errors can be estimated using Eq. III-4: 


AQ 


AQ 


n+1 


n+1 


3Q 


n+1 


3C 


3Q 


n+1 


3a 


AC = 0.297(0.025) = 0.007 


Aa = 0.183(0.05) = 0.009 


(III-20) 

v 

(III-21) 


AQ 
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AQ 


n+1 


3Q 


n+1 


3 AT 


3Q 


n+1 


3k 


A (AT) = 0.0458(0.01) = 0.0005 
Ak = -0.3866(0.15) = -0.058 


(III-22) 

(III-23) 


A comparison of the values of Eqs. III-20 to III-23 indicate that 
the largest error for this data and watershed would occur because of 
imprecision in k, with the error due to AT being insignificant. 
Obviously, the effect of error in the variables would be different for 
bther values of any of the inputs (i.e., C, a, AT, k, T, S, P) . 

While Eqs. III-20 to III-23 provide a measure of the effect of 
errors, the derivatives can be used to compute the relative sensitivity 
(Eq. III-12) . For the data of the error analysis, the relative 
sensitivities of C, a, AT, and k are: 
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(III — 24) 
(III-25) 
(III-26) 
(III-27) 


The values of Eqs. III-24 to III-27 are dimensionless and indicate 
that k and C are the most important variables, with a being only 
moderately important and AT have very low importance. 
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SECTION IV 

PROBABILISTIC MODELING 

Probabilistic modeling is an extension of sensitivity analysis 
with the following important distinctions: 1) it is not limited to 

the univariate form, 2) it allows for the interaction of the input 
variables, and 3) it provides the distribution of the design variable 
(i.e., dependent or criterion variable) and not just a single-valued 
measure of the dispersion of the design variable. Of course, there 
is a price to pay for the additional information. First, we must 
know the distributions of the input .variables , and second, the compu- 
tational requirement for a probabilistic analysis is much greater 
than that required for a linear, univariate sensitivity analysis. 

A probabilistic analysis is based on an iterative analysis of the 
design model for a sufficient number of conditions that defines the 
distribution of the design variable. While a numerical analysis of 
sensitivity requires two iterations of the design model in order to 
define the sensitivity by Eq. III-7, the probabilistic approach involves 
numerous solutions of the design model. While the sensitivity analysis 
of Eq. III-7 is univariate, with each variable being perturbed 
independently of the others, probabilistic analysis uses a simultaneous 
perturbation of all of the design variables. While the sensitivity 
analysis requires only a measure of the dispersion or error variation 
of the input variables, a probabilistic analysis requires knowledge 
of the entire distribution function of the input variables. 

Actually, probabilistic analysis is quite simple. Given the 
distribution functions of the input variables, values of the input 
variables are generated using values of random variables having the 
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distribution funct ions and used to compute a single value of the design 
variable. This process is repeated until a sufficient number of values 
of the design variables has been generated to define the distribution 
function of the design variable. 

A simple computational example can be used to illustrate the 
methodology. Let's assume a design process consists of two input 
variables, and X 2 , which are both normally distributed. The 
variable X^ has a mean and standard deviation of 5 and 1, respectively. 
Given that both have a normal distribution, values of X^ and X 2 
can be generated using a random normal number generator with the 
appropriate statistics. A sample of 25 was generated for each of the 
input variables (see Table IV- 1) . The value of a design variable Y 
is related to X.^ and X 2 by 

Y = 7.0 + 0.6XJ + 1.6X 2 (IV-1) 

The 25 values of Y generated with Eq. IV- 1 are also given in Table 
IV-1. It is Known from theory that the sum of m independent random 
variables, each normally distributed, is also normally distributed 
for a linear combination. Thus, since Eq. IV-1 is a linear equation 
and and X 2 are normally distributed, then Y must also be normally 
distributed. It can also be shown from theory that the mean and 
variance of Y are sums of the linear combinations of the input 
variables X^ and X 2 . However, the important point here is that we 
used randomly generated values of the input variables to generate the 
values of the design variable, from which an estimate of the underlying 
distribution of Y was determined. In the example used, the true 
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TABLE IV-1. Generated Values of a Design Variable (Y) as a 
Function of Two Input Variables (Xj and X.,) 


1 

5.080 

4 . 889 

17.870 

2 

1.874 

4.884 

15.939 

3 

2.574 

5.948 

IS. Obi 

4 

4.796 

4.309 

16.772 

5 

2.726 

5.952 

18.159 

6 

0.830 

3.422 

12.975 

7 

- 0.736 

5.677 

15.642 

8 

3.122 

3.740 

1 1.857 

9 

1.932 

3.270 

1 3 . 39 1 

10 

3.508 

5.421 

17.778 

11 

2.640 

>.702 

1 1.507 

12 

2.618 

5.705 

1 7 . 699 

13 

3.322 

4.956 

16.923 

14 

2.754 

4.740 

16.256 

15 

2.570 

4.842 

16 . 289 

16 

5.506 

5.321 

18.817 

17 

2.708 

5.285 

17.081 

18 

0.414 

4.444 

14.359 

19 

2.468 

5.920 

17.953 

20 

- 0.638 

5.245 

15.009 

21 

4.338 

4.239 

16.385 

22 

4.620 

5.461 

18.510 

23 

0.608 * 

3.746 

13.358 

24 

4.422 

6.241 

19.639 

25 

2.264 

2.892 

12.986 
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distribution could have been derived from theory only because the 
input variables had a normal distribution and Eq. IV- 1 defines Y to 
be a linear combination of and X In practice, the input variables 
are frequently non-normal and the relationships between Y and the 
input variables are not linear. In that case, it would be impossible 
to derive the distribution of the design variable Y from theory. 

PROBABILISTIC ANALYSIS OF SNOWMELT-RUNOFF MODEL 

The sensitivity analysis indicated that errors in the coefficients 

C, a, AT, and k can introduce significant errors into the computed 

runoff, withthe value of a, C, and k being the most important for 

the example provided. For a value of Q .of 0.503, the individual 

n+l 

errors in absolute value of 0.058, 0.009, and 0.007 for k, a, and C, 
respectively, are 11.5, 1.8,^ and 1.4 percent. These errors assume 
that each of the other coefficients are known with certainty, which 
is certainly not the case. When the errors in the coefficients 
occur simultaneously, then the accuracy of the computed snowmelt 
runoff will be considerably less than the accuracy when only one 
coefficient is in error. 

Distributions of the Coefficients 

A necessary input to a probabilistic analysis is the distribution 
functions of the variables, which would be the coefficients for the 
snowmelt runoff model. The data base was very insufficient to determine 
the distribution functions of the coefficients; therefore, the normal 
distribution was assumed. This assumption will probably not be 
critical to the conclusions, although one could make a good argument 
that for values of the recession coefficient near either zero or one, 
the error distribution of the recession coefficient would be highly 


skewed . 
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Assuming that each of the four coefficients have a normal distri- 
bution, then it is only necessary to identify the location and scale 
parameters of the distribution. Using the expected value theorem, 
the mean value provided in the User's Manual (Martinec, et al,, 1983) 
will be used to represent the location parameter. The standard 
errors derived for the sensitivity analysis will be used to represent 
the scale parameters of the distributions. These values of the 
location and scale parameters completely define the assumed normal 
populations of the four coefficients. 

Distribution of the Error in the Snowmelt Runoff 

The usual procedure for decisionmaking is to compute an average, 
or expected value, of the runoff. The predicted value represents the 
best estimate of the snowmelt runoff. While the expected value is 
important, it is also important in decisionmaking to know the likely 
variation in the predicted value due to error in the input. Such 
information can be used to compute either confidence intervals on the 
computed value or tolerance limits on the distribution of the output. 

Using the data for the ’Dischma watershed for 1974 (Martinec, 
et al., 1983), the mean one-day snowmelt runoff (m'Vs) was computed 
for each day from April 2 to July 31 and an error function was 
computed. The error about the mean predicted value was computed; this 
error distribution indicates the error due solely to error in the 
coefficients. Bias in the predicted values will not affect this error 
distribution. 

The probabilistic modeling approach was used to assess the error 
distribution. The snowmelt-runoff model of Eq. 1-1 was used with the 
Dischma data for 1974 and the one day runoff was computed for the 121 
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day period. A sequence was computed using the mean values of the 
coefficients; this sequence would be the predicted values assuming no 
error in the coefficients. Then random numbers were generated using 
a random number generator, and the values of the coefficients were 
computed using the probability distributions of the coefficients. A 
total of 250 sequences were generated, and the mean daily error was 
computed for each sequence. For the 250 mean daily errors, the mean 
of the means was 0.0115 m'Vsec and the standard deviation of the mean 
errors was 0.0468 m^/S. The distribution of the means is shown in 
Fig. IV- 1. For a 1 percent level of significance, the mean of 
0.0115 m^/S is statistically different from zero; therefore, the model 
produces slightly biased estimates for the data used. This bias is 
probably due to the high recession coefficient in the early days of the 
sequence; during this period, the distribution of the recession co- 
efficient had to be truncated in order to maintain rationality. The 
truncation has only a minor effect on the distribution of the errors, 
as shown in Fig. IV- 1. From this figure it is evident that the 
distribution is approximately normally distributed. 

Confidence Intervals on the Error of the Mean 

The above statistic can be used to derive confidence intervals on 
the mean daily estimated runoff for the 121 day period from April 2 to 
July 31. The two-sided 95 percent confidence interval would be 
±■.0918 m^/S. The two-sided 99 percent confidence interval is ±0.121 m^/S. 
These confidence intervals are valid for the mean error for the 121 
day period and not the mean error for an individual day during that 
period. 




f(e) 
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Also, these confidence intervals reflect only the error due to 
the inability to accurately predict the values of the four input 
coefficients; they do not reflect the error due to the inability of 
the model structure to represent the physical processes being modeled. 
Such errors would be substantially greater. 

The computed confidence intervals indicate that the error in the 
mean daily discharge due to the inaccuracy of the coefficients is not 
very significant. Thus, for mean long-term discharge estimates, the 
accuracy of the coefficients should not be viewed as a limiting 
factor in transferring the coefficients to adjacent watersheds in 
which the physical processes are similar. 
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SECTION V 
CONCLUSIONS 

Existing snowmelt runoff models, such as the frequently used 
regression models, are used to make single-valued estimates of runoff 
for some duration; for example, regression models have provided accurate 
estimates of snowmelt runoff volumes for periods of 30 days or more. 
While the single-value estimate can be used in decisionmaking, it 
does not provide a measure of the risk, i.e., the accuracy of the 
estimate. Inaccuracies are introduced by the structure of the model, 
the empiricism of the coefficients, the inaccuracy of the- measured 
input data, and the physical processes that are ignored in the formula- 
tion of the model. The inaccuracy of the empirical coefficients are 
of special concern because models are often calibrated on one water- 
shed and transferred to another watershed to make forecasts. Inaccuracy 
of the coefficients will result in errors in forecasts. Therefore, the 
study was conducted to show how such sources of inaccuracy can be 
evaluated. The emphasis was placed on the methodology rather than 
analysis of data since data for any one watershed are limited. 

The methodology for assessing the accuracy of forecasts consisted 
of a probabilistic modeling analysis and a sensitivity/error analysis. 
The snowmelt runoff model (Martinec, et al. , 1983) was used because it 
has a physically based structure, the input data requirements are 
more diverse, the coefficients have a physical basis, and the model has 
been applied successfully for forecasts for durations ranging from one 
day to three months or more. While the coefficients of the snowmelt 
runoff model are physically based, they must be considered as random 
variables since they cannot be determined exactly. 
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Sensitivity equations were derived for each of the four coefficients 
the runoff coefficient (C) , the degree-day factor (a) , the temperature 
adjustment (AT) , and the recession coefficient (k) . The sensitivity 
equations are given by Eqs. III-16 to III-19. The sensitivity 
equations can be used to provide a measure of the relative importance 
of the coefficients or as part of an error analysis. The results of 
a simple error analysis showed that the approximate standard error of 
the coefficients produced the largest error in the recession coefficient, 
with the errors for the runoff coefficient and the degree-day factor 
to be much less significant. The error in the runoff due to the error 
in the value of AT was very insignificant. Of course, these results 
are for one set of data. The effect of error in the coefficients is 
dependent on the watershed, the duration of the forecasts, the levels 
of the input variables (i.e., P, T, and S) , and the estimated values 
of the coefficients. However, the methodology presented here can be 
used to perform a sensitivity analysis for any case. 

The sensitivity equations show that the coefficients are inter- 
dependent and that error in one coefficient will affect the sensitivity 
of another coefficient. This suggests that a probabilistic analysis, 
should be made in addition to the sensitivity analysis. A probabilistic 
analysis was conducted for the Dischma watershed with the data for 1974. 
The results indicated that errors in the coefficients produced a 
distribution of errors in the forecasts that was approximately normally 
distributed; however, there was a slight positive bias in the mean. The 
results suggest that the input coefficients are reasonably accurate and 
that the error in the coefficients should not limit the transferability 
of the model to nearby watersheds. " 
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